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Abstract 


A numerical method is developed to study the effect of the compression deformation of the gas diffusion layer (GDL) on the performance of 
the proton exchange membrane fuel cell (PEMFC). The GDL compression deformation, caused by the clamping force, plays an important role 
in controlling the performance of PEMFC since the compression deformation affects the contact resistance, the GDL porosity distribution, and 
the cross-section area of the gas channel. In the present paper, finite element method (FEM) is used to first analyze the ohmic contact resistance 
between the bipolar plate and the GDL, the GDL deformation, and the GDL porosity distribution. Then, finite volume method is used to analyze 
the transport of the reactants and products. We investigate the effects of the GDL compression deformation, the ohmic contact resistivity, the air 
relative humidity, and the thickness of the catalyst layer (CL) on the performance of the PEMFC. The numerical results show that the fuel cell 
performance decreases with increasing compression deformation if the contact resistance is negligible, but an optimal compression deformation 


exists if the contact resistance is considerable. 
© 2007 Elsevier B.V. All rights reserved. 
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1. Introduction 


The proton exchange membrane fuel cell (PEMFC) system 
is considered as an alternative power source for electric vehicles 
by virtue of its high-energy efficiency, pollution-free character- 
istics, as well as the simplicity in its design and operation. It 
has received more and more attention in the world during the 
last decade. The membrane electrode assembly (MEA), con- 
sisting of a proton exchange membrane, a catalyst layer, and a 
gas diffusion layer, is one of the most important parts in a fuel 
cell system. Desirable functions of the GDL include effectively 
transporting the reactant gas into the catalyst layer, removing 
the liquid water from the catalyst layer to the gas channel, and 
conducting electrons with low resistance. The GDL compres- 
sion deformation caused by the clamping force influences the 
performance of the PEMFC through changing the permeability 
and conductivity. The GDL material under the bipolar plate land 
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gives a larger compression deformation than that under the chan- 
nel. This inhomogeneous compression deformation may cause 
significant changes in the local physical properties of the GDL, 
such as a nonuniform porosity distribution. The interfacial con- 
tact resistance between the GDL and the bipolar plate, controlled 
by the clamping force, is also one of the important parameters 
affecting the operational voltage loss in the PEMFC. However, 
when the interfacial contact resistivity is very small, the contact 
resistance between the GDL and the land can be neglected [1]. 
Experimental observations show that the clamping force 
plays an important role in controlling the performance of 
PEMFC [2-4]. Lee et al. [2] measured the performances of 
PEMFCs with a serpentine flow field for three types of the 
GDLs (ELAT®, TORAY™, a combination of TORAY™ and 
CARBEL®) and found that the optimal clamping force depends 
on the GDL types. Ge et al. [3] studied the effect of the GDL 
compression deformation on the PEMFC performance. They 
designed a special test fixture of fuel cells to directly measure 
the compression ratio. Nitta et al. [4] concluded that the inho- 
mogeneous compression of the GDL may lead to a significant 
local variation of the mass and charge transport properties. 
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Numerical studies have been rapidly increased in order to 
improve the performance and cost-effectiveness of PEMFCs. 
Various numerical models have been developed to investigate 
the effect of the design and operation parameters on the PEMFC 
performance [5]. In our previous work [1,6], we considered the 
inhomogeneous compression deformation of the GDL, and stud- 
ied the effects of the clamping force on the contact resistance and 
the performance of the fuel cell with an interdigitated flow field. 
The objective of our present work is to study the effects of GDL 
compression deformation on the performance of PEMFC with a 
conventional (parallel or serpentine flow pattern) flow field. We 
first use a finite element method (FEM) to analyze the contact 
pressure, the elastic deformation of the MEA, and porosity dis- 
tribution. We then analyze the gas transport using a two-phase 
unsaturated flow model. 

The numerical results show that a high compression ratio 
results in a lower porosity in the GDL, which will increase the 
transport resistance of both the gas and liquid phases. When the 
compression ratio is increased from 0% to 50%, the liquid water 
saturation in the GDL is increased by 32% and the oxygen con- 
centration decreases by 70% under the land. On the other hand, 
a high compression ratio decreases the interfacial contact resis- 
tance, and hence minimizes the electrical resistance loss inside 
the fuel cell. These two effects may work together to yield an 
optimal compression ratio, depending on the properties of the 
gas diffusion layers and the bipolar plates. Sometimes, a negli- 
gible contact resistance occurs between the bipolar plate and the 
GDL. In such a case the fuel cell performance decreases with 
the GDL compression deformation. The influence of the GDL 
compression deformation on the fuel cell performance depends 
on the operation conditions, the GDL property and the thickness 
of the catalyst layer. 


2. Model description 


The computational domain (the region surrounded by the 
dashed line as shown in Fig. la) consists of the gas diffusion 
layer, the catalyst layer, the membrane, and a period of the 
conventional gas distributor. 


2.1. FEM model 


The FEM analysis model, as shown in Fig. 1b, is used to 
analyze the GDL elastic deformation and the porosity distribu- 
tion. The deformations of the catalyst layer and the membrane 
are neglected due to their very larger elastic modulus than the 
GDL [7]. Normal displacements at the left (B6) and the right 
(B7) of the model, and the normal displacement at interface B3 
are restricted. The bipolar plate is treated as a rigid structure. 
The vertical displacement of the land is given downwards, sim- 
ulating the stacking process. The problem mentioned above can 
be described as a large deformation contact problem with the 
plane-strain simplicity. The plane-strain elements are used in 
the structure analysis and the contact elements are built at the 
possible contact region, as shown in Fig. 1b. 

The GDL material is a carbon paper, a type of non-woven 
paper of carbon fibers, which can be treated as an orthotropic 
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Fig. 1. Schematic diagrams of (a) PEMFC electrode and (b) FEM analysis 
model. 


porous material. The parameters employed in the present paper 
are given in Table 1. The GDL elastic modulus is based on the 
measurements of the stress-strain relationships given by Mishra 
et al. [8] and Ihonen [9]. 


2.2. Porosity distribution and interfacial resistance 
calculation 


The calculation methods of the porosity distribution, the GDL 
deformation, and the interfacial resistance have been described 
in details in our previous work [1]. Based on the large deforma- 
tion theory of solid mechanics, the GDL porosity after stacking 
is 


p- vo _ odo = WU- ole) a 
V Voe®v e&v 

Table 1 

GDL parameters used in the FEM analysis 

Initial porosity 0.78 


Through-plane: (10 MPa), 
in-plane: (1 GPa), References: 
[8,9] 

Case A: p =0.648(p/E,)~°7, 
Case B: p=21.74(p/E.)-°, 
Case C: p= 13.33(p/E,)!" 


Young’s modulus, E 


Interfacial resistivity p (mQ cm?) 
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where Vo is the initial total volume of the porous media, ¢o the 
initial porosity, £y the true bulk strain (logarithmic bulk strain), 
and V, is the pore volume. 

The total interfacial contact resistance in a period of the gas 


distributor is [1] 
i -1 
Aj 
($4) o 


where Rj, A;, and p; are the resistance, the contact area, and the 
resistivity of the ith contact element, respectively. The factor 
2 in Eq. (2) is due to the same contact resistances in the anode 
and cathode. Three cases of the representative contact resistivity 
listed in Table 1 are considered in the present paper. These cases 
(A, B, and C) are based on the experimental observations of 
Mishra et al. [8], Cho et al. [10], and Davies et al. [11], for 
the bipolar plates made from graphite, composite material, and 
stainless steel, respectively. The interfacial contact resistivity 
increases in the order of case A<case B <case C. 


2.3. Unsaturated flow model of the cathode 


The present simulation is based on the two-dimensional two- 
phase model developed by Lin and Nguyen [12]. The basic 
equations and parameters, including the liquid water satura- 
tion, are developed in consideration of the dependence of the 
liquid water saturation on the capillary pressure. An improved 
two-dimensional spherical-agglomerate catalyst model [13] is 
used in our analysis. The following assumptions are used in the 
two-phase flow model of the cathode to study the effect of the 
clamping force on the PEMFC performance: (1) the liquid flow 
is steady, isothermal and laminar flow; (2) the liquid flow in the 
electrode is treated as a continuous phase described by Darcy’s 
law; (3) the gas transport can be described by the diffusion mech- 
anism; and (4) the ionic potential is set to zero at the interface 
of the anode catalyst layer and the membrane. Based on these 
assumptions, the gas species conservation equations are: 


V-Nj;—-Rj)=V-(—DigVCj) — Rj =0 (3) 


where N; is the flux of the gas species, R; the source or sink terms, 
C; the concentration of the gas species, j =O or j= v denoting 


oxygen or water vapor. The effective diffusion coefficient Diy 
is 


. . IT by 
Dee = D'AT) f(b) g(s) = D! (=) fogs) (4) 


where T is the local temperature, T, the reference temperature, 
h(T) the temperature dependence function [14], b4 equals 2.334 
and 1.832 for vapor water and oxygen, respectively, and Dİ is 
the diffusivity at T, when @=1 and s=0, f) and g(s) are the 
normalized functions correlated with porosity ġ and saturation 
s, respectively. Based on the work of Tomadakis and Sotirchos 
[16], Nam and Kaviany [15] gave: 


= b2 
FO) = 4 ($ - se) as (5) 
P 


p =0.11 is a percolation threshold, the power index b is 0.521 
and 0.785, and b3 is 2 and 4, for in-plane and cross-plane dif- 
fusion, respectively. The values of these parameters are based 
on the numerical analysis results [15]. For the CL, a type of 
the packed spherical particle media, Bruggemann relation is 
employed: 


f@)g(s) = A -s)'? (6) 


Ry is the transfer rate of water at the interface of liquid water 
and water vapor: 


Ry = koll — s)yw(RT) | (PS — yw P®)q 
+ kybspwM,, | (PS — yw PE) — q) (7) 
with a switch function 


1+ lyw P8 = Pe / (yw PE = pat) 
q= 
2 


where yw = Cw/C8, s is the saturation of the liquid water in the 
pore, R is the gas constant, and My is the molecular weight of 
water. When the water partial pressure (ywP®) is greater than 
the water saturation pressure (P**') the switch function g=1. 
Otherwise, q =0. 

The liquid phase conservation equation is 


(8) 


V-Ny+R=0 (9) 


where Nw is the flux of liquid water. In the CL, there is another 
liquid water source, —2 Rọ, , due to the oxygen reduction reac- 
tion (ORR). The consumption rate of oxygen within the CL, 
Ro,, is described by the following equation [13,17,18]: 


= Co, RT 1 


(Tagg + 6)8 ]7' 
Ro, = 
H Erkre(1 a py 


Gagelagg DN 


(10) 


Eq. (10) is based on the spherical-agglomerate model, where 
H is the Henry’s constant, ragg the agglomerate radius, agg the 
effective specific surface area of the agglomerate, DN the oxy- 
gen diffusivity in the electrolyte (Nafion), kre the reaction rate 
constant, and E; is the reaction effectiveness factor for spherical 
agglomerate. E; is given as [19]: 


1 1 1 
E, = (11) 
PL \tanh3O_) 3%, 


where ®; is usually known as Thiele’s modulus for chemical 
reactions: 


Fagg krc 
QL = 12 
L=7 D, (12) 


The effective oxygen diffusivity in the agglomerate, D, is 
related to the electrolyte (Nafion) fraction in the agglomerate, 
and is given by the following equation: 


eff — pN.1.5 
DË = DNelS (13) 


where €ag¢ is the electrolyte (Nafion) fraction in the agglomerate. 
The reaction rate constant is based on Butler—Volmer equation 
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accounting for catalyst microstructure effects [13]: 


aff jref 
kre = Pt 0 
4F(— @%) } | cst 


cF 1 — ao) F 
x GA TEn) exp (| -2o ")| (14) 


where 7 is the local overpotential, ae the active or effective 
platinum surface area, and F is the Faraday’s constant. 

The transport of liquid water in the GDL and CL is driven 
by the capillary pressure. Based on the uniform gas pressure 
assumption, the flux of liquid water can be described as 


= PwK vr 


N, 
> My bw 


VPS (15) 
where pw is the density of the liquid water, and uw is the vis- 
cosity. The effective diffusion coefficient Ky}, is given by 


K Se = Ku(s)K™ (16) 


where Ky(s), representing the relative permeability of liquid 
phase [20], is proportional to the cube of the liquid water satu- 
ration, i.e.: 


Kal) = 8° (17) 


K” is the absolute permeability of water. Kozeny—Carman rela- 
tion is used to calculate the absolute permeability when the 
porosity is changed. Based on the constant Kozeny—Carman 
coefficient, one has: 


3 2 
= 
Kv = xy l= oo) 
WU — $) 
where ġo and K@ are the initial porosity and the absolute per- 


meability, respectively. The capillary pressure is related to the 
liquid water saturation via a Leverette function J(s) [21]: 


(18) 


o \'2 
Pe =o cost) J(s), 
1.417s — 2.1205? + 1.26353, if @, > 90° 
1.4170. — s) — 2.12001 =a)" (19) 


+1.263(1—s)°, if A < 90° 


J(s) = 


where ø is the surface tension of liquid water and air system, 
and 6, is the contact angle of water on the GDL. 

The liquid water transport in the membrane is driven by the 
water concentration gradient and electro-osmotic drag: 


Nw = Zna — DNycN (20) 


where CN is the liquid water concentration in the membrane. 
The net electro-osmotic drag coefficient, ng, and the diffusivity 
of the liquid water, DN, in the membrane are the functions of 
the water content in the membrane: 


DN = 1076e?4!160/303-1/1)(2.563 — 0.33 + 0.02641? 
— 0.000671A7) 21) 


2.5, 250N 
22 22 Cf 
where Cf is the concentration of the fixed charge sites in the 


membrane. At interfaces B4 and B5, aN is related to the water 
vapor activity œ [12] and can be described as 


(22) 


ng = 


ifs=0, CN = C;(0.043 + 17.81æ — 39.8507 + 36.003), 


ifs > 0, Cc = 14Cf + 2.8Cfs (23) 


2.4. Charge transport 


The local overpotential 7 can be described by the following 
equation: 


n= be — 4 -U (24) 


where U is the reference electrode potential, ®e is the elec- 
tron potential, and ®, is the ion potential. The potentials of 
electron and ion can be obtained through the following charge 
conservation equation: 


V -i4@) = V + (—K4e),eff V P+(e)) = R+) (25) 


The effective conductivities of the catalyst layer for ion, K4,eff, 
and electron, Ke,eff, are described as 


1.5 
K+ eff = Td = Oo eas 


(26) 
Keet = KEL — HA1 — age)” 


where KÀ is the ionic conductivity of the membrane, and KẸ is 
the electronic conductivity of the carbon black in agglomerate. 
The source/sink terms of electron and ion existing in the CL are 
given as 


Ry = Re = 4F RO}. (27) 
2.5. Boundary condition 


There are seven boundaries or interfaces within the studied 
domain as shown in Fig. 1. The boundary conditions are listed 
in Table 2. At the interface of the GDL and the CL (B3), the 
properties of the two layers cause a discontinuity of the liquid 
saturation. The continuity assumption of the capillary pressure 
at the interface of the two layers is used for hydrophobic GDL: 


GDL gaDt a GDL CL gt si CL 
cos(@¢ ) ( or} J(s 7) = cos(65~) (Sx) J(s~~) 
(28) 
Eq. (28) is inapplicable for hydrophilic GDL. In such a case we 
assume that s} =0 at interface B3. 


2.6. Solution procedure 


The equations together with the boundary conditions listed in 
Table 2 and described in the previous sections are solved with the 
finite element method and the finite volume method. The code is 
designed using Fortran Language. The studied domain is divided 


P. Zhou, C.W. Wu / Journal of Power Sources 170 (2007) 93-100 97 


Table 2 
The boundary conditions in the numerical analysis [12] 
Interface Variables 
C; (j: O2 or vapor) sS D, Pe cN 
B1 Cj= ce s=" N/A ielĝPL = 0 N/A 
B2 N,|SPr = 0 NwlGPL = 0 N/A Pe = N/A 
B3 Nev = N;Șt NylGOL = Nyl it = 0 Ae aC /ay = 0 
B4 Nj% = 0 Nye = Ny [pem =a = Cla s) 
B5 N/A N/A N/A N/A CN = CN(qinedey 
B6 and B7 aC;/dx=0 ds/dx = 0 d®@,/dx=0 d®,/dx =0 acN /ax = 0 


into rectangle elements with four nodes. Finite element analysis 
for the elastic deformation of the GDL is first carried out to 
obtain the real shape of the GDL after deformation. The grid is 
then changed into the approximate orthogonal curvilinear grid. 
Based on this new grid, finite volume method is used to analyze 
the distributions of the reactants and products. The parameters 
used in the calculation are listed in Table 3. 


3. Results and discussions 


Generally speaking, the output voltage and power density of a 
fuel cell will decrease with increasing the compression due to the 
decrease of the GDL porosity. However, a contrary phenomenon 
is sometimes observed in the low compression region [2,3]. 
The main reason why the performance of PEMFC increases 
with the compression is that the bulk resistivity of the GDL 
and the interfacial contact resistance decrease with the com- 
pression deformation of the GDL. The total interfacial contact 
resistance includes those coming from the bipolar plate/GDL 
interface, the GDL/CL interface, and the CL/membrane inter- 
face. However, the possible contact resistance occurring at the 
GDL/CL interface and the CL/membrane interface is considered 
to be negligible when compared with the plate/GDL interface 
resistance. 

Fig. 2 gives the numerical results with the three cases of the 
interfacial resistivity as described in Table 1. When the compres- 
sion displacement of the GDL is fixed at 0.0185 mm, the contact 
resistances for cases A, B, and C are 48.8 mQ cm, 667 mQ cm 
and 1648 mQ cm, respectively. There is a considerable influence 


—m— CASE A 
—e— CASE B 


Compression: 0.037mm(10%) 
CL Thickness: 0.016mm 


GDL Contact Angle:100° 


Power Density (J cm”) 


‘0.0 0.2 0.4 0.6 0.8 1.0 
Current Density (A cm) 


Fig. 2. Polarization curves for the contact resistivity described by three different 
equations. 


of the contact resistance on the output voltage and power den- 
sity of the fuel cell. It can be inferred that the contact resistance 
mainly affects the slope of the polarization curve at the ohmic 
region. 

Fig. 3 shows the polarization curves at different compres- 
sion ratios for case C. The compression ratio is defined as the 


Table 3 
The physical properties of MEA and employed parameters in the calculations 


Parameter Value 


Gas diffusion layer properties 


Channel width (half) 0.05 cm 

Shoulder width (half) 0.05 cm 

Initial porosity (Toray 120), og 0.78 

Initial thickness (Toray 120), Ago 0.037 cm 

Absolute permeability, Kapt [22] 8.69 x 1078 cm? 

Contact angle, 6°! [23] 100°; 89° 
Electronic conductivity «GPL 

In plane (Toray 120) 213Scm7! 

Through plane (Toray 120) 12.5 S cm7! 


Catalyst layer properties 
Initial porosity, ost [13] 0.1 


Initial thickness, nce 

Absolute permeability, Kc. 

Contact angle, ot 

Effective platinum surface area, agit [13] 

Electrolyte fraction in agglomerate, agg 
[24] 

Carbon black conductivity [12] 

Henry’s constant, H [25] 

Ref. O2 concentration, cH [26] 


Membrane properties 


Initial thickness, hy 


Ionic conductivity, «"°™ 


Operation conditions 


Temperature, T 

Cathode gas pressure 

Mole fraction of oxygen in the air inlet 
Humidity in the air inlet 

Humidity at anode 

Water saturation level at channel 
Evaporation rate constant, ky 
Condensation rate constant, ke 
Cathode transfer coefficient, œe [26] 


ref 


Reference exchange current density, ip 


[26] 


Reference electrode potential, U 


0.0016 cm; 0.003 cm 
107!! cm? 

110° 

303797.5 cm? cm~? 
0.5 


100S cm7! 
0.3125 atm m? mol”! 
0.85 mol m~? 


0.012 cm 
0.1Scm7! 


65°C 

l atm 

0.21 

75%; 100% 

100% 

0 

100 atm~! s7! 

100 mol atm~! s~! cm~? 
1 (=0.8 V) 

0.61 (<0.8 V) 

3.85 x 1078 A cm~? (>0.8 V) 


1.5 x 1076 A cm™°? (<0.8 V) 
1.1 V 
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Fig. 3. Polarization curves for several compression ratios as indicated. 


ratio of the normal compression displacement to the initial GDL 
thickness. We find that the voltage increases with increasing 
compression due to the decrease in the ohmic overpotential. It 
is shown that there exists an optimal compression displacement 
(~0.0555 mm). The limiting current decreases with an increase 
in the GDL compression deformation. 

The dependence of the limiting current on the transport ability 
of the reactant gas varies with the porosity and the liquid water 
distribution. Fig. 4 shows the polarization curves with different 
cathode air humidities in order to study the influence of liquid 
water saturation. With the decrease in the relative humidity, the 
saturation of liquid water also decreases due to the interfacial 
transfer of water between the liquid water and the water vapor. 
When the relative humidity is decreased to ~75%, the liquid 
water saturation in most of the regions of the GDL deceases to 
about null. From Fig. 4, it is shown that the higher the humidity, 
the lower the limiting current and the greater the compression 
effect. 

There is a significant effect of the liquid water saturation 
on the PEMFC performance. The saturation of liquid water is 
affected by the GDL wettability. Gurau et al. [23] have studied 
the internal contact angle of water within the GDL material. It 
is estimated to be 89° for 30% PTFE and 100° for 70% PTFE. 
Comparing what is shown in Figs. 4 and 5, it is found that the 
limiting current of the fuel cell decreases from ~0.98 A cm? to 
~0.66 A cm? for 15% of the compression ratio when the contact 
angle decreases from 100° to 89°. The effect of the compression 


1.0% CL Thickness: 0.0016cm 
GDL Contact Angle:100° 


0.8 


Voltage (V) 
oO 
O 


RH Compression 
—# 100% 0.0555mm(15%) 
—a— 100% 0.1850mm(50%) 
-#@ - 75% 0.0555mm(15%) 
-@- 75% 0.1850mm(50%) 


>o 
a 
Power Density (J cm’) 


0.0 0.2 0.4 0.6 0.8 1.0 1.2) 
Current Density (A cm’) 


Fig. 4. Polarization curves with different cathode air humidities. 
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Fig. 5. Polarization curves for different compressions and different contact 
resistivity models: (a) contact resistivity described by case A and (b) contact 
resistivity described by case C. 


on the fuel cell performance is greater when the contact angle 
is 89° than when the contact angle is 100°. Comparing what is 
shown in Figs. 5a and b, the limiting current density decreases 
with increasing compression. However, there exists an optimal 
compression ratio to give the highest power density when case 
C is used. But no apparent optimal compression ratio is found 
for case A. Generally speaking, there always exists an optimal 
compression ratio for a practical FEMFC system because the 
interfacial contact resistivity approaches infinite when the com- 
pression ratio is close to zero (zero interfacial pressure). For 
case A, the optimal compression ratio may range from 0 to 5%. 
The experimental measurements of Ge et al. [3] show that the 
fuel cell performance decreases with an increase in the com- 
pression ratio when the compression ratio ranges from 15% to 
45%. However, it is estimated that the optimal compression ratio 
may be between 0% and 15%, which was not be found in their 
experiments. 

The GDL compression always decreases the porosity, thus 
reducing the transport ability of the reactant gas and liquid water 
in the GDL. Fig. 6 shows the distribution of the mole fraction 
of oxygen within the GDL and CL and the distribution of the 
liquid water saturation within the GDL when the compression 
of the GDL is set at 0% and 50%. A higher water saturation may 
exist in the GDL when the compression ratio is 50% due to the 
reduced porosity. From Fig. 6b and Fig. 6d, the maximum liquid 
water saturation is 0.489 and 0.647 when the compression ratio 
is set at 0% and 50%, respectively. It is difficult to transport 
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(a) Land 


Fig. 6. Distributions of (a) oxygen mole fraction with 0% compression ratio, 
(b) liquid water saturation level with 0% compression ratio, (c) oxygen mole 
fraction with 50% compression ratio and (d) liquid water saturation level with 
50% compression ratio (cell voltage: 0.3 V; contact angle of the GDL: 89°; 
contact resistance: neglected). 


oxygen from the channel/GDL interface to the CL due to the 
reduced effective void for the gas diffusion. When the compres- 
sion ratio increases from 0% to 50%, the oxygen concentration 
decreases by about 70% under the land. In this paper, we use 
a given boundary condition at the channel/GDL interface for 
the liquid water saturation and the gas components. However, 
the GDL compression affects the cross-section area of the chan- 
nel which affects the removal of the liquid water and the gas 
near the channel/GDL interface. When the compression ratio 
reaches 25% and 50%, the GDL surface goes into the chan- 
nel by about 0.0552 mm and 0.1115 mm, respectively, due to its 
elastic deformation. The influence of the cross-section area of 
the channel on the performance of the PEMFC can be studied 
in a three-dimensional flow field analysis. 

Finally, we study the combined effects of the CL thickness 
and GDL contact angle on the performance of the fuel cell as 
shown in Fig. 7a and Fig. 7b. Usually, a larger CL thickness 
gives a higher limiting current and a more significant effect of 
the compression on the performance, especially when the GDL 
contact angle is 89°. 
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Fig. 7. Polarization curves at different compressions and CL thickness: (a) the 
contact angle is 89° and (b) the contact angle is 100°. 


What is shown above can explain both of the experimental 
observations given by Lee et al. [2] and Ge et al. [3]. However, 
the effect of the compression on the fuel cell performance in the 
present simulation is slightly smaller than what was measured 
by Ge et al. [3]. One possible reason for this may be that the 
effects of compression on the physical and geometrical proper- 
ties of the membrane and catalyst layers were not considered in 
the present simulation. Furthermore, the GDL compression will 
also affect the anode flow field. This may cause some difference 
between the experiments and the numerical simulation. It should 
be stressed that the contact resistivity depends on not only the 
land conduction property but also the GDL conduction property. 
That is why Lee et al. [2] obtained different results for different 
GDLs in their experiments. 


4. Conclusions 


In order to study the effect of the GDL compression defor- 
mation on the PEMFC performance, a simulation method 
is developed in the present paper. Finite element method 
(FEM) and two-phase flow model are used to investigate 
the polarization curves of the PEMFC when the GDL defor- 
mation in the cathode with the conventional flow field is 
considered. We have studied the combined effects of the 
compression ratio, the interfacial contact resistivity, the CL 
thickness, the GDL wettability, and the air relative humidity. 
From the numerical results, we have reached the following 
conclusions: 
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(1) The GDL compression deformation decreases the contact 
resistance, which in turn causes a drop in the ohmic overpo- 
tential of the polarization curve, especially for case B and 
case C (high interfacial contact resistivity). Until now, most 
of the experimental observations for the contact resistivity 
are higher than what is described by case B. Therefore, the 
contact resistance plays an important role in polarization 
curve simulation. 

(2) Large GDL compression deformation and large porosity 
variation reduce the transport ability of the reactant gas 
and liquid water in the GDL. On the other hand, the GDL 
compression deformation reduces the ohmic overpotential. 
The optimal compression ratio results from the combined 
effects of the GDL compression deformation on the contact 
resistance and the mass transfer. 

(3) The effect of the GDL compression deformation is con- 
siderable large in the case of high current density, high air 
humidity, and thick catalyst layer due to the increasing liquid 
water saturation. 
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